correlation_part2
Comparing newt and snake phenotype correlation
newt and snake data from Hanifin et al. 2008
snewt <- read.csv("~/Desktop/snake_newt_real_data.csv", header = TRUE)
snewt$mean_res = snewt$mean_toxin+snewt$dif
cor(snewt$mean_toxin, snewt$mean_res)## [1] 0.7660207(sites <- data.frame(longitude = snewt$lon, latitude = snewt$lat))## longitude latitude
## 1 124.4964 49.9415
## 2 124.5653 49.7453
## 3 125.4494 49.6506
## 4 120.4164 37.2972
## 5 124.2026 41.7558
## 6 121.9018 37.8534
## 7 121.5566 37.0030
## 8 121.7326 40.5401
## 9 120.5733 38.5813
## 10 124.0598 41.2864
## 11 122.3255 37.5630
## 12 121.1893 35.6440
## 13 120.5724 34.7420
## 14 123.3556 39.4096
## 15 123.6314 40.9396
## 16 116.6764 46.7360
## 17 123.3874 44.6282
## 18 124.3873 42.7068
## 19 122.4554 42.1040
## 20 121.8897 44.2164
## 21 123.5642 43.0962
## 22 123.9238 46.1651
## 23 123.9454 46.2729
## 24 124.0168 48.0405
## 25 122.6448 47.1003
## 26 123.8839 46.9770
## 27 122.3490 48.3387
## 28 122.2156 45.6487just_loc <- data.frame(snewt$lon, snewt$lat)
snewt$lon<- snewt$lon *-1
locations_sf <- st_as_sf(snewt, coords = c("lon", "lat"), crs = 4326)
mapview(locations_sf)library("rnaturalearth")
library("rnaturalearthdata")
library("tools")
library("maps")##
## Attaching package: 'maps'## The following object is masked from 'package:purrr':
##
## mapworld <- ne_countries(scale = "medium", returnclass = "sf")
states <- st_as_sf(maps::map("state", plot = FALSE, fill = TRUE))
states$ID <- toTitleCase(states$ID)
ggplot(data = world) +
geom_sf() +
geom_sf(data = states, fill = NA)+
geom_point(data = snewt, aes(x = lon, y = lat, fill = mean_res), size = 6, shape = 22, ) +
geom_point(data = snewt, aes(x = lon, y = lat, fill = mean_toxin), size = 4, shape = 21, ) +
scale_fill_viridis_c()+theme(axis.text.x = element_text(angle = 45, hjust=1, size = 8))+
coord_sf(xlim = c(-127, -115), ylim = c(34, 50.5), expand = FALSE) 
ggplot(data = world) +
geom_sf() +
geom_sf(data = states, fill = NA)+
#geom_point(data = snewt, aes(x = lon, y = lat, fill = mean_res), size = 6, shape = 22, ) +
geom_point(data = snewt, aes(x = lon, y = lat, fill = mean_toxin), size = 4, shape = 21, ) +
scale_fill_viridis_c()+theme(axis.text.x = element_text(angle = 45, hjust=1, size = 8))+
coord_sf(xlim = c(-127, -115), ylim = c(34, 50.5), expand = FALSE)
ggplot(data = world) +
geom_sf() +
geom_sf(data = states, fill = NA)+
geom_point(data = snewt, aes(x = lon, y = lat, fill = mean_res), size = 6, shape = 22, ) +
#geom_point(data = snewt, aes(x = lon, y = lat, fill = mean_toxin), size = 4, shape = 21, ) +
scale_fill_viridis_c()+theme(axis.text.x = element_text(angle = 45, hjust=1, size = 8))+
coord_sf(xlim = c(-127, -115), ylim = c(34, 50.5), expand = FALSE)
More data from the slim simulation
might be a specific generation or a small set of generations, If I can get the data to be saved and then interpreted by R make a list in a column how to proportion and recreate what I see in the slim map (the different box and their dims)
GA1 experiment values:
- 1e-8*(0.005)^2 = 2.5e-13
- 1e-9*(0.05)^2 = 2.5e-12
- 1e-10*(0.5)^2 = 2.5e-11
- 1e-11*(5.0)^2 = 2.5e-10
GA2 experiment values:
- 1e-8*(0.05)^2 = 2.5e-11
- 1e-9*(0.1581139)^2 = 2.5e-11
- 1e-10*(0.5)^2 = 2.5e-11
- 1e-11*(1.581139)^2 = 2.5e-11
- 1e-12*(5)^2 = 2.5e-11
#cor_grid_30 <- list.files(path="~/Desktop/GA1_rep", pattern="_grid_", full.names=TRUE, recursive=FALSE)
#cor_grid_30_G2 <- list.files(path="~/Desktop/GA2_rep", pattern="_grid_", full.names=TRUE, recursive=FALSE)
cor_grid_30 <- list.files(path="~/Desktop/GA1_csv_notperarea", pattern="_grid_", full.names=TRUE, recursive=FALSE)
cor_grid_30_G2 <- list.files(path="~/Desktop/GA2_csv_notperarea", pattern="_grid_", full.names=TRUE, recursive=FALSE)
cor_grid_30_file<-read.csv(cor_grid_30[1], header = TRUE)
replica <- unlist(str_split(unlist(str_split(cor_grid_30[1],"rep_"))[3], "_"))[1]
cor_grid_30_file$rep<-replica
cor_grid_30_file$sim <-"A"
cor_grid_30_file$rep <-unlist(str_split(unlist(str_split(cor_grid_30[1],"rep_"))[3], "_"))[1]
cor_grid_30_file$snake_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30[1],"snake_mu_rate_"))[2], "_"))[1]
cor_grid_30_file$newt_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30[1],"newt_mu_rate_"))[2], "_"))[1]
cor_grid_30_file$snake_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30[1],"snake_mu_effect_sd_"))[2], "_"))[1]
cor_grid_30_file$newt_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30[1],"newt_mu_effect_sd_"))[2], "_"))[1]
cor_grid_30_G2_file<-read.csv(cor_grid_30_G2[1], header = TRUE)
replica <- unlist(str_split(unlist(str_split(cor_grid_30_G2[1],"rep_"))[3], "_"))[1]
cor_grid_30_G2_file$rep<-replica
cor_grid_30_G2_file$sim <-"A"
cor_grid_30_G2_file$rep <-unlist(str_split(unlist(str_split(cor_grid_30_G2[1],"rep_"))[3], "_"))[1]
cor_grid_30_G2_file$snake_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30_G2[1],"snake_mu_rate_"))[2], "_"))[1]
cor_grid_30_G2_file$newt_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30_G2[1],"newt_mu_rate_"))[2], "_"))[1]
cor_grid_30_G2_file$snake_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30_G2[1],"snake_mu_effect_sd_"))[2], "_"))[1]
cor_grid_30_G2_file$newt_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30_G2[1],"newt_mu_effect_sd_"))[2], "_"))[1]
oldstand=0
letter=1
for(i in seq(1, length(cor_grid_30_G2), 1)){
stand <- unlist(str_split(unlist(str_split(cor_grid_30_G2[i],"newt_mu_rate_"))[2], "_"))[1]
if (stand!=oldstand){
oldstand=stand
#letter = letter + 1
}
letter = i
GA2_lit_temp_file <- read.csv(cor_grid_30_G2[i], header = TRUE)
GA2_lit_temp_file$sim <-LETTERS[letter]
GA2_lit_temp_file$rep <-unlist(str_split(unlist(str_split(cor_grid_30_G2[i],"rep_"))[3], "_"))[1]
GA2_lit_temp_file$snake_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30_G2[i],"snake_mu_rate_"))[2], "_"))[1]
GA2_lit_temp_file$newt_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30[i],"newt_mu_rate_"))[2], "_"))[1]
GA2_lit_temp_file$snake_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30[i],"snake_mu_effect_sd_"))[2], "_"))[1]
GA2_lit_temp_file$newt_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30[i],"newt_mu_effect_sd_"))[2], "_"))[1]
cor_grid_30_G2_file <- rbind(cor_grid_30_G2_file, GA2_lit_temp_file)
if (i<=length(cor_grid_30)){
cor_temp_file <- read.csv(cor_grid_30[i], header = TRUE)
cor_temp_file$sim <-LETTERS[letter]
cor_temp_file$rep <-unlist(str_split(unlist(str_split(cor_grid_30[i],"rep_"))[3], "_"))[1]
cor_temp_file$snake_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30[i],"snake_mu_rate_"))[2], "_"))[1]
cor_temp_file$newt_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30[i],"newt_mu_rate_"))[2], "_"))[1]
cor_temp_file$snake_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30[i],"snake_mu_effect_sd_"))[2], "_"))[1]
cor_temp_file$newt_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30[i],"newt_mu_effect_sd_"))[2], "_"))[1]
cor_grid_30_file <- rbind(cor_grid_30_file, cor_temp_file)
rep <-unlist(str_split(unlist(str_split(cor_grid_30[i],"rep_"))[3], "_"))[1]
snake_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30[i],"snake_mu_rate_"))[2], "_"))[1]
newt_mu_rate <- unlist(str_split(unlist(str_split(cor_grid_30[i],"newt_mu_rate_"))[2], "_"))[1]
snake_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30[i],"snake_mu_effect_sd_"))[2], "_"))[1]
newt_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_grid_30[i],"newt_mu_effect_sd_"))[2], "_"))[1]
print(paste0("Simulation ", LETTERS[letter],": Snake mu-rate & effect sd (", snake_mu_rate,", ", snake_mu_effect_sd, ") Newt mu-rate & effect sd (", newt_mu_rate, ", ", newt_mu_effect_sd, ")", " rep: ", rep))
}
}## [1] "Simulation A: Snake mu-rate & effect sd (1.0e-08, 0.005) Newt mu-rate & effect sd (1.0e-08, 0.005) rep: 0"
## [1] "Simulation B: Snake mu-rate & effect sd (1.0e-08, 0.005) Newt mu-rate & effect sd (1.0e-09, 0.05) rep: 0"
## [1] "Simulation C: Snake mu-rate & effect sd (1.0e-08, 0.005) Newt mu-rate & effect sd (1.0e-10, 0.5) rep: 0"
## [1] "Simulation D: Snake mu-rate & effect sd (1.0e-08, 0.005) Newt mu-rate & effect sd (1.0e-11, 5.0) rep: 0"
## [1] "Simulation E: Snake mu-rate & effect sd (1.0e-09, 0.05) Newt mu-rate & effect sd (1.0e-08, 0.005) rep: 0"
## [1] "Simulation F: Snake mu-rate & effect sd (1.0e-09, 0.05) Newt mu-rate & effect sd (1.0e-09, 0.05) rep: 0"
## [1] "Simulation G: Snake mu-rate & effect sd (1.0e-09, 0.05) Newt mu-rate & effect sd (1.0e-10, 0.5) rep: 0"
## [1] "Simulation H: Snake mu-rate & effect sd (1.0e-09, 0.05) Newt mu-rate & effect sd (1.0e-11, 5.0) rep: 0"
## [1] "Simulation I: Snake mu-rate & effect sd (1.0e-10, 0.5) Newt mu-rate & effect sd (1.0e-08, 0.005) rep: 0"
## [1] "Simulation J: Snake mu-rate & effect sd (1.0e-10, 0.5) Newt mu-rate & effect sd (1.0e-09, 0.05) rep: 0"
## [1] "Simulation K: Snake mu-rate & effect sd (1.0e-10, 0.5) Newt mu-rate & effect sd (1.0e-10, 0.5) rep: 0"
## [1] "Simulation L: Snake mu-rate & effect sd (1.0e-10, 0.5) Newt mu-rate & effect sd (1.0e-11, 5.0) rep: 0"
## [1] "Simulation M: Snake mu-rate & effect sd (1.0e-11, 5.0) Newt mu-rate & effect sd (1.0e-08, 0.005) rep: 0"
## [1] "Simulation N: Snake mu-rate & effect sd (1.0e-11, 5.0) Newt mu-rate & effect sd (1.0e-09, 0.05) rep: 0"
## [1] "Simulation O: Snake mu-rate & effect sd (1.0e-11, 5.0) Newt mu-rate & effect sd (1.0e-10, 0.5) rep: 0"
## [1] "Simulation P: Snake mu-rate & effect sd (1.0e-11, 5.0) Newt mu-rate & effect sd (1.0e-11, 5.0) rep: 0"cor_grid_30_file$snake <- paste0(cor_grid_30_file$snake_mu_rate, "_",cor_grid_30_file$snake_mu_effect_sd)
cor_grid_30_file$newt <- paste0(cor_grid_30_file$newt_mu_rate, "_",cor_grid_30_file$newt_mu_effect_sd)
cor_grid_30_file$newt_snake_rep <- paste0("Newt_",cor_grid_30_file$newt_mu_rate, "_", cor_grid_30_file$newt_mu_effect_sd, "_Snake_",cor_grid_30_file$snake_mu_rate, "_", cor_grid_30_file$snake_mu_effect_sd, "_", cor_grid_30_file$rep)
#GA1_mean = cor_grid_30_file[ , grepl( "mean" , names( cor_grid_30_file ) ) ]
#unique(cor_grid_30_file$newt_snake_rep)
GA1<-cor_grid_30_file[which(cor_grid_30_file$generation>=10000 & cor_grid_30_file$gen<=110000 ),]
#GA2<-cor_grid_30_G2_file[which(cor_grid_30_G2_file$generation>=10000 & cor_grid_30_G2_file$gen<=110000 ),]
#length(GA1$generation)
#GA1<-cor_grid_30_file[which(cor_grid_30_file$generation<=50 & cor_grid_30_file$newt_snake_rep=="Newt_1.0e-08_0.005_Snake_1.0e-09_0.05_2"),]
GA1_mean = GA1[,grepl("generation|mean", names(GA1)) ]
#How to remake the dataframe into gen, newt_value, snake_value, location (x, y)?
GA1_t <- data.frame(t(GA1_mean))
loc_x = 1
loc_y = 5
skip = 1
GA1_t$loc_x = 0
GA1_t$loc_y = 0
GA1_t$species = 0
grid_info <- as.character(rownames(GA1_t))
for (i in 2:nrow(GA1_t)){
#print(as.numeric(unlist(str_split(i,"_pheno_"))[2]))
GA1_t$species[i] = unlist(str_split(grid_info[i],"_mean"))[1]
if(loc_y==0){
loc_x = loc_x + 1
loc_y = 5
}
#print(paste0(loc_x, " ", loc_y))
GA1_t$loc_x[i]= loc_x
GA1_t$loc_y[i] = loc_y
if (skip %% 2 == 0 ){
loc_y = loc_y - 1
}
skip = skip +1
}
GA1_t$names<- rownames(GA1_t)
GA1_mean = cor_grid_30_file[,grepl( "mean" , names(cor_grid_30_file)) ]
GA1_mean$generation <- cor_grid_30_file$generation
GA1_mean$rep <- cor_grid_30_file$rep
GA1_mean$snake_mu_rate <- cor_grid_30_file$snake_mu_rate
GA1_mean$newt_mu_rate <- cor_grid_30_file$newt_mu_rate
GA1_mean$snake_mu_effect_sd <- cor_grid_30_file$snake_mu_effect_sd
GA1_mean$newt_mu_effect_sd <- cor_grid_30_file$newt_mu_effect_sd
GA1_mean$snake <- cor_grid_30_file$snake
GA1_mean$newt <- cor_grid_30_file$newt
GA1_mean$newt_snake_rep <- cor_grid_30_file$newt_snake_rep
data_long_mean <- gather(GA1_mean, condition, measurement, newt_mean_pheno_0:snake_mean_pheno_24, factor_key=TRUE)
data_long_mean$lox <- GA1_t[match(data_long_mean$condition, GA1_t$names),"loc_x"]
data_long_mean$loy <- GA1_t[match(data_long_mean$condition, GA1_t$names),"loc_y"]
data_long_mean$species <- GA1_t[match(data_long_mean$condition, GA1_t$names),"species"]
GA1_max = cor_grid_30_file[,grepl( "max" , names(cor_grid_30_file)) ]
GA1_max$generation <- cor_grid_30_file$generation
GA1_max$rep <- cor_grid_30_file$rep
GA1_max$snake_mu_rate <- cor_grid_30_file$snake_mu_rate
GA1_max$newt_mu_rate <- cor_grid_30_file$newt_mu_rate
GA1_max$snake_mu_effect_sd <- cor_grid_30_file$snake_mu_effect_sd
GA1_max$newt_mu_effect_sd <- cor_grid_30_file$newt_mu_effect_sd
GA1_max$snake <- cor_grid_30_file$snake
GA1_max$newt <- cor_grid_30_file$newt
GA1_max$newt_snake_rep <- cor_grid_30_file$newt_snake_rep
data_long_max <- gather(GA1_max, condition, measurement, newt_max_pheno_0:snake_max_pheno_24, factor_key=TRUE)
data_long_max$lox <- data_long_mean$lox
data_long_max$loy <- data_long_mean$loy
data_long_max$species <- data_long_mean$speciesGA1_onefile<-data_long_mean[which(data_long_mean$newt_snake_rep=="Newt_1.0e-10_0.5_Snake_1.0e-10_0.5_0" & data_long_max$gen==10981 ),]
#GA1_onefile<-data_long_mean[which(data_long_mean$newt_snake_rep==unique(data_long_mean$newt_snake_rep)[16] & data_long_max$gen==29981 ),]
#GA1_onefile<-data_long_max[which(data_long_max$newt_snake_rep=="Newt_1.0e-10_0.5_Snake_1.0e-10_0.5_0" & data_long_max$gen==29051 ),]
ggplot(GA1_onefile)+
geom_point(data = GA1_onefile[which(GA1_onefile$species=="snake"),], aes(x = lox, y = loy, fill = measurement), size = 10, shape = 22, ) +
geom_point(data = GA1_onefile[which(GA1_onefile$species=="newt"),], aes(x = lox, y = loy, fill = measurement), size = 7, shape = 21, ) +
scale_fill_viridis_c()+
theme_bw()
unique(data_long_mean$newt_snake_rep)## [1] "Newt_1.0e-08_0.005_Snake_1.0e-08_0.005_0"
## [2] "Newt_1.0e-09_0.05_Snake_1.0e-08_0.005_0"
## [3] "Newt_1.0e-10_0.5_Snake_1.0e-08_0.005_0"
## [4] "Newt_1.0e-11_5.0_Snake_1.0e-08_0.005_0"
## [5] "Newt_1.0e-08_0.005_Snake_1.0e-09_0.05_0"
## [6] "Newt_1.0e-09_0.05_Snake_1.0e-09_0.05_0"
## [7] "Newt_1.0e-10_0.5_Snake_1.0e-09_0.05_0"
## [8] "Newt_1.0e-11_5.0_Snake_1.0e-09_0.05_0"
## [9] "Newt_1.0e-08_0.005_Snake_1.0e-10_0.5_0"
## [10] "Newt_1.0e-09_0.05_Snake_1.0e-10_0.5_0"
## [11] "Newt_1.0e-10_0.5_Snake_1.0e-10_0.5_0"
## [12] "Newt_1.0e-11_5.0_Snake_1.0e-10_0.5_0"
## [13] "Newt_1.0e-08_0.005_Snake_1.0e-11_5.0_0"
## [14] "Newt_1.0e-09_0.05_Snake_1.0e-11_5.0_0"
## [15] "Newt_1.0e-10_0.5_Snake_1.0e-11_5.0_0"
## [16] "Newt_1.0e-11_5.0_Snake_1.0e-11_5.0_0"GA1_onefile_2<-data_long_mean[which(data_long_mean$newt_snake_rep=="Newt_1.0e-08_0.005_Snake_1.0e-08_0.005_0" & data_long_max$gen==10981 ),]
ggplot(GA1_onefile_2)+
geom_point(data = GA1_onefile_2[which(GA1_onefile_2$species=="snake"),], aes(x = lox, y = loy, fill = measurement), size = 10, shape = 22, ) +
geom_point(data = GA1_onefile_2[which(GA1_onefile_2$species=="newt"),], aes(x = lox, y = loy, fill = measurement), size = 7, shape = 21, ) +
scale_fill_viridis_c()+
theme_bw()
snake_info<-data_long_mean[which(data_long_mean$species=="snake"),]
newt_info<-data_long_mean[which(data_long_mean$species=="newt"),]
names(newt_info)[names(newt_info) == 'measurement'] <- 'newt_toxin'
newt_info$snake_res<-snake_info$measurement
newt_info$dif<-newt_info$snake_res-newt_info$newt_toxin
pattern<-"Newt_1.0e-08_0.005_Snake_1.0e-08_0.005_0"
scatter_pheno<-newt_info[which(newt_info$newt_snake_rep==pattern & newt_info$generation>14000 & newt_info$generation<15000 ),]
#scatter_pheno<-newt_info[which(newt_info$newt_snake_rep==pattern & newt_info$generation==20981 ),]
scatter_pheno[scatter_pheno == 0] <- NA
scatter_pheno[scatter_pheno == 1] <- NA
ggplot(scatter_pheno)+
ggtitle(paste0(" pattern ", pattern, " From ", min(scatter_pheno$generation), " gen to ", max(scatter_pheno$generation), " gen"))+
geom_point(aes(y = newt_toxin, x = snake_res, color=generation), size = 3) ## Warning: Removed 350 rows containing missing values (`geom_point()`).
theme_bw()## List of 97
## $ line :List of 6
## ..$ colour : chr "black"
## ..$ linewidth : num 0.5
## ..$ linetype : num 1
## ..$ lineend : chr "butt"
## ..$ arrow : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_line" "element"
## $ rect :List of 5
## ..$ fill : chr "white"
## ..$ colour : chr "black"
## ..$ linewidth : num 0.5
## ..$ linetype : num 1
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ text :List of 11
## ..$ family : chr ""
## ..$ face : chr "plain"
## ..$ colour : chr "black"
## ..$ size : num 11
## ..$ hjust : num 0.5
## ..$ vjust : num 0.5
## ..$ angle : num 0
## ..$ lineheight : num 0.9
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ title : NULL
## $ aspect.ratio : NULL
## $ axis.title : NULL
## $ axis.title.x :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 2.75points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.x.top :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 0
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 2.75points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.x.bottom : NULL
## $ axis.title.y :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : num 90
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 2.75points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.y.left : NULL
## $ axis.title.y.right :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 0
## ..$ angle : num -90
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 2.75points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : chr "grey30"
## ..$ size : 'rel' num 0.8
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 2.2points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x.top :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 0
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 2.2points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x.bottom : NULL
## $ axis.text.y :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 1
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 2.2points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.y.left : NULL
## $ axis.text.y.right :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 2.2points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.ticks :List of 6
## ..$ colour : chr "grey20"
## ..$ linewidth : NULL
## ..$ linetype : NULL
## ..$ lineend : NULL
## ..$ arrow : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_line" "element"
## $ axis.ticks.x : NULL
## $ axis.ticks.x.top : NULL
## $ axis.ticks.x.bottom : NULL
## $ axis.ticks.y : NULL
## $ axis.ticks.y.left : NULL
## $ axis.ticks.y.right : NULL
## $ axis.ticks.length : 'simpleUnit' num 2.75points
## ..- attr(*, "unit")= int 8
## $ axis.ticks.length.x : NULL
## $ axis.ticks.length.x.top : NULL
## $ axis.ticks.length.x.bottom: NULL
## $ axis.ticks.length.y : NULL
## $ axis.ticks.length.y.left : NULL
## $ axis.ticks.length.y.right : NULL
## $ axis.line : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ axis.line.x : NULL
## $ axis.line.x.top : NULL
## $ axis.line.x.bottom : NULL
## $ axis.line.y : NULL
## $ axis.line.y.left : NULL
## $ axis.line.y.right : NULL
## $ legend.background :List of 5
## ..$ fill : NULL
## ..$ colour : logi NA
## ..$ linewidth : NULL
## ..$ linetype : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ legend.margin : 'margin' num [1:4] 5.5points 5.5points 5.5points 5.5points
## ..- attr(*, "unit")= int 8
## $ legend.spacing : 'simpleUnit' num 11points
## ..- attr(*, "unit")= int 8
## $ legend.spacing.x : NULL
## $ legend.spacing.y : NULL
## $ legend.key :List of 5
## ..$ fill : chr "white"
## ..$ colour : logi NA
## ..$ linewidth : NULL
## ..$ linetype : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ legend.key.size : 'simpleUnit' num 1.2lines
## ..- attr(*, "unit")= int 3
## $ legend.key.height : NULL
## $ legend.key.width : NULL
## $ legend.text :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : 'rel' num 0.8
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ legend.text.align : NULL
## $ legend.title :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ legend.title.align : NULL
## $ legend.position : chr "right"
## $ legend.direction : NULL
## $ legend.justification : chr "center"
## $ legend.box : NULL
## $ legend.box.just : NULL
## $ legend.box.margin : 'margin' num [1:4] 0cm 0cm 0cm 0cm
## ..- attr(*, "unit")= int 1
## $ legend.box.background : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ legend.box.spacing : 'simpleUnit' num 11points
## ..- attr(*, "unit")= int 8
## $ panel.background :List of 5
## ..$ fill : chr "white"
## ..$ colour : logi NA
## ..$ linewidth : NULL
## ..$ linetype : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ panel.border :List of 5
## ..$ fill : logi NA
## ..$ colour : chr "grey20"
## ..$ linewidth : NULL
## ..$ linetype : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ panel.spacing : 'simpleUnit' num 5.5points
## ..- attr(*, "unit")= int 8
## $ panel.spacing.x : NULL
## $ panel.spacing.y : NULL
## $ panel.grid :List of 6
## ..$ colour : chr "grey92"
## ..$ linewidth : NULL
## ..$ linetype : NULL
## ..$ lineend : NULL
## ..$ arrow : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_line" "element"
## $ panel.grid.major : NULL
## $ panel.grid.minor :List of 6
## ..$ colour : NULL
## ..$ linewidth : 'rel' num 0.5
## ..$ linetype : NULL
## ..$ lineend : NULL
## ..$ arrow : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_line" "element"
## $ panel.grid.major.x : NULL
## $ panel.grid.major.y : NULL
## $ panel.grid.minor.x : NULL
## $ panel.grid.minor.y : NULL
## $ panel.ontop : logi FALSE
## $ plot.background :List of 5
## ..$ fill : NULL
## ..$ colour : chr "white"
## ..$ linewidth : NULL
## ..$ linetype : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ plot.title :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : 'rel' num 1.2
## ..$ hjust : num 0
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 5.5points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ plot.title.position : chr "panel"
## $ plot.subtitle :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 5.5points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ plot.caption :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : 'rel' num 0.8
## ..$ hjust : num 1
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 5.5points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ plot.caption.position : chr "panel"
## $ plot.tag :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : 'rel' num 1.2
## ..$ hjust : num 0.5
## ..$ vjust : num 0.5
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ plot.tag.position : chr "topleft"
## $ plot.margin : 'margin' num [1:4] 5.5points 5.5points 5.5points 5.5points
## ..- attr(*, "unit")= int 8
## $ strip.background :List of 5
## ..$ fill : chr "grey85"
## ..$ colour : chr "grey20"
## ..$ linewidth : NULL
## ..$ linetype : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ strip.background.x : NULL
## $ strip.background.y : NULL
## $ strip.clip : chr "inherit"
## $ strip.placement : chr "inside"
## $ strip.text :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : chr "grey10"
## ..$ size : 'rel' num 0.8
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 4.4points 4.4points 4.4points 4.4points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ strip.text.x : NULL
## $ strip.text.x.bottom : NULL
## $ strip.text.x.top : NULL
## $ strip.text.y :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : num -90
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ strip.text.y.left :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : num 90
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ strip.text.y.right : NULL
## $ strip.switch.pad.grid : 'simpleUnit' num 2.75points
## ..- attr(*, "unit")= int 8
## $ strip.switch.pad.wrap : 'simpleUnit' num 2.75points
## ..- attr(*, "unit")= int 8
## - attr(*, "class")= chr [1:2] "theme" "gg"
## - attr(*, "complete")= logi TRUE
## - attr(*, "validate")= logi TRUEcor_files <- list.files(path="~/Desktop/csv_GA1_long", pattern="_cor_", full.names=TRUE, recursive=FALSE)
cor_GA2_files <- list.files(path="~/Desktop/GA2_rep", pattern="_cor_", full.names=TRUE, recursive=FALSE)
lit_30 <- list.files(path="~/Desktop/csv_GA1_long", pattern="_lit_", full.names=TRUE, recursive=FALSE)
lit_30_G2 <- list.files(path="~/Desktop/GA2_rep", pattern="_lit_", full.names=TRUE, recursive=FALSE)
lit_30_file<- data.frame()
lit_30_G2_file <- data.frame()
number=1
#cor_file <- read.csv(cor_files[number], header = TRUE)
#cor_file$sim <-"A"
#cor_GA2_file <- read.csv(cor_GA2_files[number], header = TRUE)
#cor_GA2_file$sim <-"A"
cor_file<- data.frame()
cor_GA2_file<- data.frame()
oldstand=0
letter=0
newt_unique=0
newt_numb = 0
for(i in seq(1, length(cor_GA2_files), 1)){
stand <- unlist(str_split(unlist(str_split(cor_GA2_files[i],"newt_mu_rate_"))[2], "_"))[1]
up <- unlist(str_split(unlist(str_split(cor_GA2_files[i],"newt_mu_effect_sd_"))[2], "_"))[1]
phrase = paste0(stand, "_", up)
if (stand!=oldstand){
oldstand=stand
letter = letter + 1
}
if (phrase=="1.0e-08_0.05" ){
newt_numb = 1
}
if (phrase=="1.0e-09_0.158114"){
newt_numb = 2
}
if (phrase=="1.0e-10_0.5"){
newt_numb = 3
}
if (phrase=="1.0e-11_1.58114"){
newt_numb = 4
}
if (phrase=="1.0e-12_5.0"){
newt_numb = 5
}
GA2_cor_temp_file <- read.csv(cor_GA2_files[i], header = TRUE)
GA2_cor_temp_file$sim <-LETTERS[letter]
GA2_cor_temp_file$size <-newt_numb
GA2_cor_temp_file$rep <-unlist(str_split(unlist(str_split(cor_GA2_files[i],"rep_"))[3], "_"))[1]
GA2_cor_temp_file$snake_mu_rate <- unlist(str_split(unlist(str_split(cor_GA2_files[i],"snake_mu_rate_"))[2], "_"))[1]
GA2_cor_temp_file$newt_mu_rate <- unlist(str_split(unlist(str_split(cor_GA2_files[i],"newt_mu_rate_"))[2], "_"))[1]
GA2_cor_temp_file$snake_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_GA2_files[i],"snake_mu_effect_sd_"))[2], "_"))[1]
GA2_cor_temp_file$newt_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_GA2_files[i],"newt_mu_effect_sd_"))[2], "_"))[1]
cor_GA2_file <- rbind(cor_GA2_file, GA2_cor_temp_file)
GA2_lit_temp_file <- read.csv(lit_30_G2[i], header = TRUE)
GA2_lit_temp_file$sim <-LETTERS[letter]
GA2_lit_temp_file$size <-newt_numb
GA2_lit_temp_file$rep <-unlist(str_split(unlist(str_split(lit_30_G2[i],"rep_"))[3], "_"))[1]
GA2_lit_temp_file$snake_mu_rate <- unlist(str_split(unlist(str_split(lit_30_G2[i],"snake_mu_rate_"))[2], "_"))[1]
GA2_lit_temp_file$newt_mu_rate <- unlist(str_split(unlist(str_split(lit_30_G2[i],"newt_mu_rate_"))[2], "_"))[1]
GA2_lit_temp_file$snake_mu_effect_sd <- unlist(str_split(unlist(str_split(lit_30_G2[i],"snake_mu_effect_sd_"))[2], "_"))[1]
GA2_lit_temp_file$newt_mu_effect_sd <- unlist(str_split(unlist(str_split(lit_30_G2[i],"newt_mu_effect_sd_"))[2], "_"))[1]
lit_30_G2_file <- rbind(lit_30_G2_file, GA2_lit_temp_file)
if (i<=length(cor_files)){
stand <- unlist(str_split(unlist(str_split(cor_files[i],"newt_mu_rate_"))[2], "_"))[1]
up <- unlist(str_split(unlist(str_split(cor_files[i],"newt_mu_effect_sd_"))[2], "_"))[1]
phrase = paste0(stand, "_", up)
if (phrase=="1.0e-08_0.005"){
newt_numb = 1
}
if (phrase=="1.0e-09_0.05"){
newt_numb = 2
}
if (phrase=="1.0e-10_0.5"){
newt_numb = 3
}
if (phrase=="1.0e-11_5.0"){
newt_numb = 4
}
cor_temp_file <- read.csv(cor_files[i], header = TRUE)
cor_temp_file$sim <-LETTERS[letter]
cor_temp_file$size <-newt_numb
cor_temp_file$rep <-unlist(str_split(unlist(str_split(cor_files[i],"rep_"))[3], "_"))[1]
cor_temp_file$snake_mu_rate <- unlist(str_split(unlist(str_split(cor_files[i],"snake_mu_rate_"))[2], "_"))[1]
cor_temp_file$newt_mu_rate <- unlist(str_split(unlist(str_split(cor_files[i],"newt_mu_rate_"))[2], "_"))[1]
cor_temp_file$snake_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_files[i],"snake_mu_effect_sd_"))[2], "_"))[1]
cor_temp_file$newt_mu_effect_sd <- unlist(str_split(unlist(str_split(cor_files[i],"newt_mu_effect_sd_"))[2], "_"))[1]
cor_file <- rbind(cor_file, cor_temp_file)
GA1_temp_file <- read.csv(lit_30[i], header = TRUE)
GA1_temp_file$sim <-LETTERS[letter]
GA1_temp_file$size <-newt_numb
GA1_temp_file$rep <-unlist(str_split(unlist(str_split(lit_30[i],"rep_"))[3], "_"))[1]
GA1_temp_file$snake_mu_rate <- unlist(str_split(unlist(str_split(lit_30[i],"snake_mu_rate_"))[2], "_"))[1]
GA1_temp_file$newt_mu_rate <- unlist(str_split(unlist(str_split(lit_30[i],"newt_mu_rate_"))[2], "_"))[1]
GA1_temp_file$snake_mu_effect_sd <- unlist(str_split(unlist(str_split(lit_30[i],"snake_mu_effect_sd_"))[2], "_"))[1]
GA1_temp_file$newt_mu_effect_sd <- unlist(str_split(unlist(str_split(lit_30[i],"newt_mu_effect_sd_"))[2], "_"))[1]
lit_30_file <- rbind(lit_30_file, GA1_temp_file)
#print(paste0("Simulation ", LETTERS[letter],": Snake mu-rate & effect sd (", cor_temp_file$snake_mu_rate[i],", ", cor_temp_file$snake_mu_effect_sd[i], ") Newt mu-rate & effect sd (", cor_temp_file$newt_mu_rate[i], ", ",cor_temp_file$ newt_mu_effect_sd[i], ")" ))
}
}
lit_30_file$dif <- lit_30_file$Snake_mean_Pheno-lit_30_file$Newt_mean_Pheno
lit_30_G2_file$dif <- lit_30_G2_file$Snake_mean_Pheno-lit_30_G2_file$Newt_mean_Phenoboxplot_cors_by_gen <- function(file, title, subt, x_label, y_label, min_gen, max_gen){
real_snewt_cor = cor(snewt$mean_toxin, snewt$mean_res)
file$newt_snake_rep <- paste0(file$newt_mu_rate, "_", file$newt_mu_effect_sd, "_", file$snake_mu_rate, "_", file$snake_mu_effect_sd, "_", file$rep)
file$newt_snake <- paste0(file$newt_mu_rate, "_", file$newt_mu_effect_sd, "_", file$snake_mu_rate, "_", file$snake_mu_effect_sd)
unique(file$sim)
plot <- ggplot(file[which(file$gen>=min_gen & file$gen<=max_gen),]) +
ylim(-1.1, 1.1)+
ggtitle(title, subtitle = paste0(subt, " From ", min_gen, " gen to ", max_gen, " gen"))+
xlab(x_label) + ylab(y_label)+
geom_boxplot(aes(x=factor(rep), y=as.numeric(mean_newt_pheno_By_mean_snake_pheno), color = "CA"))+
geom_hline(yintercept=real_snewt_cor, linetype="dashed")+
geom_hline(yintercept=0, linetype="solid")+
#geom_boxplot(aes(x=factor("CB"), y=as.numeric(num_newts_By_num_snakes), color = "CB"))+
#geom_boxplot(aes(x=factor("CC"), y=as.numeric(sum_newt_pheno_By_num_snake), color = "CC"))+
#geom_boxplot(aes(x=factor("CD"), y=as.numeric(sum_snake_pheno_By_num_newt), color = "CD"))+
theme(plot.title = element_text(size=5)) +
facet_wrap(~newt_snake)+
theme_bw() + theme(axis.text.x = element_text(angle = 45, hjust=1, size = 8), legend.position="none",
strip.background = element_blank(), strip.text.x = element_blank(),
title =element_text(size=8, face='bold'), axis.title.x = element_blank(),
axis.title.y = element_blank())
return(plot)
}
var_cor <- function(file, cfile, title, subt, x_label, y_label, min_gen, max_gen){
real_snewt_cor = cor(snewt$mean_toxin, snewt$mean_res)
real_snewt_var = var(snewt$dif)
cfile$newt_snake_rep <- paste0("Newt_", cfile$newt_mu_rate, "_", cfile$newt_mu_effect_sd, "_Snake_", cfile$snake_mu_rate, "_", cfile$snake_mu_effect_sd, "_", cfile$rep)
var_cor_table=data.frame()
for(i in 1:length(unique(file$newt_snake_rep))){
pattern<-unique(file$newt_snake_rep)[i]
temp_file<-file[which(file$gen>=min_gen & file$gen<=max_gen & file$newt_snake_rep==pattern),]
temp_cor_file<-cfile[which(cfile$gen>=min_gen & cfile$gen<=max_gen & cfile$newt_snake_rep==pattern),]
var_cor_table_temp<- data.frame("pattern"=pattern, "var"= var(temp_file$dif),"cor" = cor(temp_file$newt_toxin, temp_file$snake_res), "mean_corr" = mean(temp_cor_file$mean_newt_pheno_By_mean_snake_pheno))
var_cor_table <- rbind(var_cor_table, var_cor_table_temp)
}
var_cor_table_temp<- data.frame("pattern"="real", "var"= real_snewt_var, "cor" = real_snewt_cor, "mean_corr" = real_snewt_cor)
var_cor_table <- rbind(var_cor_table, var_cor_table_temp)
plot <- ggplot(var_cor_table) +
ggtitle(title, subtitle = paste0(subt, " From ", min_gen, " gen to ", max_gen, " gen"))+
xlab(x_label) + ylab(y_label)+
geom_point(aes(x=as.numeric(var), y=as.numeric(cor), color = pattern), shape=17, size=3)+
geom_point(aes(x=as.numeric(var), y=as.numeric(mean_corr), color = pattern), size=3)+
theme(plot.title = element_text(size=5)) +
theme_bw() + theme(axis.text.x = element_text(angle = 45, hjust=1, size = 8))
return(plot)
}
trends_by_gen <- function(file, cfile, title, subt, x_label, x_label_cor, y_label, min_gen, max_gen, pattern){
file$newt_snake_rep <- paste0(file$newt_mu_rate, "_", file$newt_mu_effect_sd, "_", file$snake_mu_rate, "_", file$snake_mu_effect_sd, "_", file$rep)
cfile$newt_snake_rep <- paste0(cfile$newt_mu_rate, "_", cfile$newt_mu_effect_sd, "_", cfile$snake_mu_rate, "_", cfile$snake_mu_effect_sd, "_", cfile$rep)
plot <- ggplot(file[which(file$gen>=min_gen & file$gen<=max_gen & file$newt_snake_rep==pattern),]) +
#geom_line(aes(generation, dif, group=rep, col="dif"))+
ggtitle(title, subtitle = paste0(subt," pattern ", pattern, " From ", min_gen, " gen to ", max_gen, " gen"))+
geom_line(aes(generation, Snake_mean_Pheno, group=rep, color="Snake_mean_Pheno"),)+
geom_line(aes(generation, Newt_mean_Pheno, group=rep, color="Newt_mean_Pheno"))+
#ylim(0,10)+
#xlab(x_label) + ylab(y_label)+
scale_color_manual(name = "Type",
values = c( "Snake_mean_Pheno" = "steelblue", "Newt_mean_Pheno" = "darkred"),
labels = c("Snake", "Newt"))+
theme_bw() +theme(legend.position="none",
strip.background = element_blank(), strip.text.x = element_blank(),
title =element_text(size=8, face='bold'), axis.title.x = element_blank(),
axis.title.y = element_blank(), axis.ticks.x = element_blank(), axis.text.x = element_blank(),
plot.margin=unit(c(1,1,+0.3,1), "cm"))
cor_plot <- ggplot(cfile[which(cfile$gen>=min_gen & cfile$gen<=max_gen & cfile$newt_snake_rep==pattern),]) +
geom_line(aes(generation, mean_newt_pheno_By_mean_snake_pheno, col="CA"))+
#ggtitle(title, subtitle = subt) +
ylim(-1.1, 1.1)+
#xlab(x_label_cor) + ylab(y_label)+
theme_bw() +theme(axis.text.x = element_text(angle = 45, hjust=1, size = 8), legend.position="none",
strip.background = element_blank(), strip.text.x = element_blank(),
title =element_text(size=8, face='bold'), axis.title.x = element_blank(),
axis.title.y = element_blank(),
plot.margin=unit(c(-0.5,1,1,1), "cm"))
file<-file[which(file$gen>=min_gen & file$gen<=max_gen & file$newt_snake_rep==pattern),]
cfile<-cfile[which(cfile$gen>=min_gen & cfile$gen<=max_gen & cfile$newt_snake_rep==pattern),]
toDelete <- seq(0, nrow(cfile), 2)
cfile_short <- cfile[-toDelete, ]
figure <- ggarrange(plot, cor_plot + font("x.text", size = 10),heights = c(2, 1.5),ncol = 1, nrow = 2, align = "v")
print(paste0("pattern ", pattern))
print(paste0("Cor between average snake pheno and local cor ", cor(file$Snake_mean_Pheno, cfile_short$mean_newt_pheno_By_mean_snake_pheno)))
print(paste0("Cor between average newt pheno and local cor ", cor(file$Newt_mean_Pheno, cfile_short$mean_newt_pheno_By_mean_snake_pheno)))
print(paste0("Cor between average dif pheno and local cor ", cor(file$dif, cfile_short$mean_newt_pheno_By_mean_snake_pheno)))
return(figure)
}
four_types_boxplot_cors_by_gen <- function(file,title, subt, x_label, y_label, min_gen, max_gen){
file$newt_snake <- paste0(file$newt_mu_rate, "_", file$newt_mu_effect_sd, "_", file$snake_mu_rate, "_", file$snake_mu_effect_sd)
plot <- ggplot(file[which(file$gen>=min_gen & file$gen<=max_gen),]) +
ylim(-1.1, 1.1)+
ggtitle(title, subtitle = paste0(subt, " From ", min_gen, " gen to ", max_gen, " gen"))+
xlab(x_label) + ylab(y_label)+
geom_boxplot(aes(x=factor("CA"), y=as.numeric(mean_newt_pheno_By_mean_snake_pheno), color = "CA"))+
geom_boxplot(aes(x=factor("CB"), y=as.numeric(num_newts_By_num_snakes), color = "CB"))+
geom_boxplot(aes(x=factor("CC"), y=as.numeric(sum_newt_pheno_By_num_snake), color = "CC"))+
geom_boxplot(aes(x=factor("CD"), y=as.numeric(sum_snake_pheno_By_num_newt), color = "CD"))+
geom_hline(yintercept=0, linetype="solid")+
theme(plot.title = element_text(size=5)) +
facet_wrap(~newt_snake)+
theme_bw() + theme(axis.text.x = element_text(angle = 45, hjust=1, size = 8), legend.position="none",
strip.background = element_blank(), strip.text.x = element_blank(),
title =element_text(size=8, face='bold'), axis.title.x = element_blank(),
axis.title.y = element_blank())
return(plot)
}boxplot_cors_by_gen(cor_file, "Simulation Correlation Comparied to Real Correlation", "GA1", "newt/snake GAs&rep", "Cor", 4000, 5000)
four_types_boxplot_cors_by_gen(cor_file, "Simulation four Correlation", "GA1", "newt/snake GAs&rep", "Cor", 4000, 5000)
boxplot_cors_by_gen(cor_GA2_file, "Simulation Correlation Comparied to Real Correlation", "GA2", "newt/snake GAs&rep", "Cor", 4000, 5000)
trends_by_gen(lit_30_file, cor_file, "test", "GA1", "newt/snake GAs&rep","other", "Cor", min_gen=0, max_gen=5000, pattern="1.0e-10_0.5_1.0e-10_0.5_0")## [1] "pattern 1.0e-10_0.5_1.0e-10_0.5_0"
## [1] "Cor between average snake pheno and local cor 0.257236703085727"
## [1] "Cor between average newt pheno and local cor 0.0460065772310917"
## [1] "Cor between average dif pheno and local cor 0.0735100676282853"
trends_by_gen(lit_30_file, cor_file, "test", "GA1", "newt/snake GAs&rep","other", "Cor", min_gen=10000, max_gen=15000, pattern="1.0e-10_0.5_1.0e-10_0.5_0")## [1] "pattern 1.0e-10_0.5_1.0e-10_0.5_0"
## [1] "Cor between average snake pheno and local cor -0.0244914135160876"
## [1] "Cor between average newt pheno and local cor 0.331733410464774"
## [1] "Cor between average dif pheno and local cor -0.208485161185655"
trends_by_gen(lit_30_file, cor_file, "test", "GA1", "newt/snake GAs&rep","other", "Cor", min_gen=0, max_gen=5000, pattern="1.0e-08_0.005_1.0e-08_0.005_0")## [1] "pattern 1.0e-08_0.005_1.0e-08_0.005_0"
## [1] "Cor between average snake pheno and local cor 0.103828794302331"
## [1] "Cor between average newt pheno and local cor 0.299141985390285"
## [1] "Cor between average dif pheno and local cor -0.372704622259776"
trends_by_gen(lit_30_file, cor_file, "test", "GA1", "newt/snake GAs&rep","other", "Cor", min_gen=0, max_gen=15000, pattern="1.0e-09_0.05_1.0e-08_0.005_0")## [1] "pattern 1.0e-09_0.05_1.0e-08_0.005_0"
## [1] "Cor between average snake pheno and local cor 0.0153139143652654"
## [1] "Cor between average newt pheno and local cor -0.0284561836126627"
## [1] "Cor between average dif pheno and local cor 0.0280010323520237"
trends_by_gen(lit_30_file, cor_file, "test", "GA1", "newt/snake GAs&rep","other", "Cor", min_gen=0, max_gen=15000, pattern="1.0e-08_0.005_1.0e-09_0.05_0")## [1] "pattern 1.0e-08_0.005_1.0e-09_0.05_0"
## [1] "Cor between average snake pheno and local cor 0.0298135650824854"
## [1] "Cor between average newt pheno and local cor 0.0977110523590969"
## [1] "Cor between average dif pheno and local cor -0.00538798771324545"
trends_by_gen(lit_30_G2_file, cor_GA2_file, "test", "GA1", "newt/snake GAs&rep","other", "Cor", min_gen=0, max_gen=15000, pattern="1.0e-08_0.05_1.0e-08_0.05_0")## [1] "pattern 1.0e-08_0.05_1.0e-08_0.05_0"
## [1] "Cor between average snake pheno and local cor 0.117586087993649"
## [1] "Cor between average newt pheno and local cor 0.0720273011182107"
## [1] "Cor between average dif pheno and local cor 0.0342810647748496"
snake_info<-data_long_mean[which(data_long_mean$species=="snake"),]
newt_info<-data_long_mean[which(data_long_mean$species=="newt"),]
names(newt_info)[names(newt_info) == 'measurement'] <- 'newt_toxin'
newt_info$snake_res<-snake_info$measurement
newt_info$dif<-newt_info$snake_res-newt_info$newt_toxin
var_cor(newt_info, cor_file, "var by cor", "mini", "var", "cor", 10000, 15000)